dafsar1<-function(n){
library(mvtnorm)

y<-array(,dim=c(6,n))
e<-array(,dim=c(6,n))

f<-array(,dim=c(2,n))
v<-array(,dim=c(2,n))

#parameter:
f0<-array(c(0,0),dim=c(2,1))

h<-array(c(1,1,1,0,0,0,
           0,0,0,1,1,1),dim=c(6,2))

fp<-array(c(.8,0,0,.8),dim=c(2,2))
r<-array(c(1,0,0,0,0,0,
           0,1,0,0,0,0,
           0,0,1,0,0,0,
           0,0,0,1,0,0,
           0,0,0,0,1,0,
           0,0,0,0,0,1),dim=c(6,6))
qq<-array(c(.36,.18,.18,.36),dim=c(2,2))

for (i in 1:n){
   v[,i]=t(rmvnorm(1,c(0,0),qq))
   f[,i]=fp%*%f0+v[,i]
   e[,i]=t(rmvnorm(1,c(0,0,0,0,0,0),r))
   y[,i]=h%*%f[,i]+e[,i]
   f0<-f[,i]
}
return(list(y=t(y)[1001:n,],f=t(f)[1001:n,],v=t(v)[1001:n,],e=t(e)[1001:n,]))
}